This dataset contains many different animals across space and time, recording a variety of data that seems meaningful to the groups. In this project I chose to work with primates.
The primates in this dataset contain Adult Body Mass, and Adult Head-to-Body lengths, which is what I have chose to model, along with using their taxonomic groups. The hopes are to model the relationship of the animal’s size(lengths) to their mass, while considering their groups in a multilevel model.
library(brms)
## Loading required package: Rcpp
## Loading 'brms' package (version 2.12.0). Useful instructions
## can be found by typing help('brms'). A more detailed introduction
## to the package is available through vignette('brms_overview').
##
## Attaching package: 'brms'
## The following object is masked from 'package:stats':
##
## ar
library(dagitty)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(reshape2)
library(ggplot2)
library(RColorBrewer)
#Load Dataset
df.animals <- read.csv2("data/BoneLengths.csv", sep=",", fileEncoding="UTF-8-BOM",as.is=T) #Encoding fixes name ERR
df.animals
#What order has the most data?
#counts <- df.primate %>% group_by(df.primate$MSW05_Order) %>% count() %>% arrange(desc(n))
counts <- df.animals %>% group_by(df.animals$MSW05_Genus) %>% count() %>% arrange(desc(n))
counts
df.animals$X5.1_AdultBodyMass_g <- as.numeric(df.animals$X5.1_AdultBodyMass_g)
df.animals$X8.1_AdultForearmLen_mm <- as.numeric(df.animals$X8.1_AdultForearmLen_mm)
df.animals$X13.1_AdultHeadBodyLen_mm <- as.numeric(df.animals$X13.1_AdultHeadBodyLen_mm)
Primates were chosen from the dataset
df.animals <- df.animals %>% filter(df.animals$X5.1_AdultBodyMass_g != -999)
df.animals <- df.animals %>% filter (df.animals$X5.1_AdultBodyMass_g != 1)
#df.animals <- df.animals %>% filter (df.animals$X8.1_AdultForearmLen_mm != 1)
df.animals <- df.animals %>% filter (df.animals$X13.1_AdultHeadBodyLen_mm != 1)
df.animals <- df.animals %>% filter (df.animals$X13.1_AdultHeadBodyLen_mm != -999)
df.primate <- df.animals %>% filter(df.animals$MSW05_Order == "Primates")
#df.bats <- df.animals %>% filter(df.animals$MSW05_Order)
#Filter by Reference
#df.primate <- df.primate %>% mutate(found_ref = grepl("2655", References))
#df.primate <- df.primate %>% filter(found_ref == TRUE)
counts <- df.primate %>% group_by(df.primate$MSW05_Genus) %>% count() %>% arrange(desc(n))
#print(counts)
df.primate
#Dags
primate.dag <- dagitty("dag {
Order -> Family
Family -> Genus
Genus -> Species
Order -> Forearm_Length
Order -> Head_Body_Length
Order -> Body_Mass
Family -> Forearm_Length
Family -> Head_Body_Length
Family -> Body_Mass
Genus -> Forearm_Length
Genus -> Head_Body_Length
Genus -> Body_Mass
Species -> Forearm_Length
Species -> Head_Body_Length
Species -> Body_Mass
Forearm_Length -> Body_Mass
Head_Body_Length -> Body_Mass
Forearm_Length <-> Head_Body_Length
}")
plot(graphLayout (primate.dag))
adjustmentSets(primate.dag, exposure = c("Order","Family", "Genus", "Species", "Head_Body_Length"), outcome = "Body_Mass", effect="total")
impliedConditionalIndependencies(primate.dag)
## Family _||_ Species | Genus
## Genus _||_ Order | Family
## Order _||_ Species | Genus
## Order _||_ Species | Family
This dag was narrowed down due to the lack of data in each species and no measurments of forearm length.
primate.dag <- dagitty("dag {
Order -> Family
Family -> Genus
Order -> Head_Body_Length
Order -> Body_Mass
Family -> Head_Body_Length
Family -> Body_Mass
Genus -> Head_Body_Length
Genus -> Body_Mass
Head_Body_Length -> Body_Mass
}")
plot(graphLayout (primate.dag))
adjustmentSets(primate.dag, exposure = c("Order","Family", "Genus", "Species", "Head_Body_Length"), outcome = "Body_Mass", effect="total")
## {}
impliedConditionalIndependencies(primate.dag)
## Genus _||_ Order | Family
#df.primate <- df.primate %>% mutate(scaled_AdultForearmLen_mm = scale(df.primate$X8.1_AdultForearmLen_mm))
df.primate <- df.primate %>% mutate(scaled_Adult_Mass_g = scale(df.primate$X5.1_AdultBodyMass_g))
df.primate <- df.primate %>% mutate(scaled_Adult_Length_g = scale(df.primate$X13.1_AdultHeadBodyLen_mm))
#Visualize the data
ggplot(df.primate, aes(y=X5.1_AdultBodyMass_g, x=X13.1_AdultHeadBodyLen_mm , color=MSW05_Family)) + geom_point() + coord_cartesian(xlim = c(100,900), ylim = c(-5,50000))
ggplot(df.primate, aes(y=X5.1_AdultBodyMass_g, x=X13.1_AdultHeadBodyLen_mm , color=MSW05_Family)) + geom_point()
ggplot(df.primate, aes(y=X5.1_AdultBodyMass_g, x=X13.1_AdultHeadBodyLen_mm , color=MSW05_Genus)) + geom_point() + coord_cartesian(xlim = c(100,900), ylim = c(-5,50000)) +theme(legend.position="none")
ggplot(df.primate, aes(y=X5.1_AdultBodyMass_g, x=X13.1_AdultHeadBodyLen_mm , color=MSW05_Genus)) + geom_point() + theme(legend.position="none")
This visualization shows that there are clusters for each Primate Family as well as each Genus. These will be used to model a hierachical relationship with the head-body lengths in a multilevel regression.
#Priors and Model
get_prior(scaled_Adult_Mass_g ~ X13.1_AdultHeadBodyLen_mm + (1+X13.1_AdultHeadBodyLen_mm|MSW05_Family/MSW05_Genus) , data=df.primate, family = "gaussian")
my.priors = prior("normal(0,1)", class= "b" ) + prior("normal(0,1)", class = "Intercept") + prior("normal(0,1)", class="sd")
fit.prior <- brm(scaled_Adult_Mass_g ~ scaled_Adult_Length_g + (1+scaled_Adult_Length_g|MSW05_Family/MSW05_Genus), data=df.primate, family = "gaussian", iter=2000, chains = 4, cores = 4, prior = my.priors, sample_prior = "only" )
## Compiling the C++ model
## Start sampling
plot(fit.prior, ask=F)
pp_check(fit.prior)
## Using 10 posterior samples for ppc type 'dens_overlay' by default.
summary(fit.prior)
## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: scaled_Adult_Mass_g ~ scaled_Adult_Length_g + (1 + scaled_Adult_Length_g | MSW05_Family/MSW05_Genus)
## Data: df.primate (Number of observations: 179)
## Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup samples = 4000
##
## Group-Level Effects:
## ~MSW05_Family (Number of levels: 15)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.80 0.60 0.03 2.17 1.00
## sd(scaled_Adult_Length_g) 0.80 0.61 0.03 2.24 1.00
## cor(Intercept,scaled_Adult_Length_g) -0.00 0.56 -0.94 0.93 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 3818 2435
## sd(scaled_Adult_Length_g) 4060 2381
## cor(Intercept,scaled_Adult_Length_g) 8141 2279
##
## ~MSW05_Family:MSW05_Genus (Number of levels: 67)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.81 0.62 0.03 2.28 1.00
## sd(scaled_Adult_Length_g) 0.81 0.59 0.04 2.23 1.00
## cor(Intercept,scaled_Adult_Length_g) -0.00 0.58 -0.96 0.95 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 3608 1752
## sd(scaled_Adult_Length_g) 3849 2157
## cor(Intercept,scaled_Adult_Length_g) 9277 1898
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept -0.03 0.97 -1.94 1.81 1.00 9305
## scaled_Adult_Length_g -0.01 0.99 -1.97 1.94 1.00 10641
## Tail_ESS
## Intercept 2839
## scaled_Adult_Length_g 2748
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 10.88 12.66 0.33 40.61 1.00 5840 1989
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
This prior model starts out assuming there is no relationship between head-body length and mass.
plot(marginal_effects(fit.prior), ask=F)
## Warning: Method 'marginal_effects' is deprecated. Please use
## 'conditional_effects' instead.
#Vespertilionidae Eptesicus nasutus
df.pred <- data.frame(scaled_Adult_Length_g = seq(-2,2,0.1), MSW05_Family="Cercopithecidae", MSW05_Genus="Cercocebus" )
df.pred
pred <- predict(fit.prior, newdata = df.pred)
df.pred.all <-cbind(df.pred, pred)
ggplot(df.pred.all, aes(x=scaled_Adult_Length_g,y=Estimate)) + geom_point()
fit.bones <- brm(scaled_Adult_Mass_g ~ scaled_Adult_Length_g + (1+scaled_Adult_Length_g|MSW05_Family/MSW05_Genus), data=df.primate, family = "gaussian", iter=2000, chains = 4, cores = 4, prior = my.priors)
## Compiling the C++ model
## recompiling to avoid crashing R session
## Start sampling
## Warning: There were 1 divergent transitions after warmup. Increasing adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning: Examine the pairs() plot to diagnose sampling problems
summary(fit.bones)
## Warning: There were 1 divergent transitions after warmup. Increasing adapt_delta
## above 0.8 may help. See http://mc-stan.org/misc/warnings.html#divergent-
## transitions-after-warmup
## Family: gaussian
## Links: mu = identity; sigma = identity
## Formula: scaled_Adult_Mass_g ~ scaled_Adult_Length_g + (1 + scaled_Adult_Length_g | MSW05_Family/MSW05_Genus)
## Data: df.primate (Number of observations: 179)
## Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
## total post-warmup samples = 4000
##
## Group-Level Effects:
## ~MSW05_Family (Number of levels: 15)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.37 0.10 0.21 0.58 1.00
## sd(scaled_Adult_Length_g) 0.38 0.08 0.25 0.57 1.00
## cor(Intercept,scaled_Adult_Length_g) 0.98 0.04 0.87 1.00 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1226 1726
## sd(scaled_Adult_Length_g) 1591 1979
## cor(Intercept,scaled_Adult_Length_g) 1253 1339
##
## ~MSW05_Family:MSW05_Genus (Number of levels: 67)
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sd(Intercept) 0.10 0.02 0.06 0.14 1.00
## sd(scaled_Adult_Length_g) 0.09 0.02 0.06 0.14 1.00
## cor(Intercept,scaled_Adult_Length_g) 0.91 0.09 0.64 1.00 1.00
## Bulk_ESS Tail_ESS
## sd(Intercept) 1558 2329
## sd(scaled_Adult_Length_g) 1595 2321
## cor(Intercept,scaled_Adult_Length_g) 1494 1833
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
## Intercept -0.17 0.11 -0.39 0.04 1.00 814
## scaled_Adult_Length_g 0.33 0.11 0.11 0.55 1.00 817
## Tail_ESS
## Intercept 1087
## scaled_Adult_Length_g 1259
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.11 0.01 0.10 0.12 1.00 3349 2930
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
plot(fit.bones)
pp_check(fit.bones)
## Using 10 posterior samples for ppc type 'dens_overlay' by default.
pp_check(fit.bones, type="intervals")
## Using all posterior samples for ppc type 'intervals' by default.
This assumes the total regression, so we explore each regression by grouped clusters.
plot(marginal_effects(fit.bones), ask=F, points=TRUE)
## Warning: Method 'marginal_effects' is deprecated. Please use
## 'conditional_effects' instead.
Here, we see each multilevel regression by family. As we know, there
coef(fit.bones)$MSW05_Family
## , , Intercept
##
## Estimate Est.Error Q2.5 Q97.5
## Aotidae -0.37695315 0.21457876 -0.7966475 0.05137255
## Atelidae -0.02754708 0.05499033 -0.1322643 0.08533033
## Cebidae -0.28592281 0.10098433 -0.4848308 -0.09027649
## Cercopithecidae -0.04878720 0.03132087 -0.1080289 0.01668930
## Cheirogaleidae -0.34761766 0.16432183 -0.6873702 -0.03250054
## Daubentoniidae -0.22211031 0.18518926 -0.5861527 0.13844094
## Galagidae -0.36537515 0.15421333 -0.6869757 -0.07050006
## Hominidae 0.91745711 0.20267376 0.4860554 1.30213604
## Hylobatidae -0.14171765 0.06170918 -0.2577011 -0.01282146
## Indriidae -0.19026778 0.07725071 -0.3481673 -0.03949486
## Lemuridae -0.31928781 0.07109031 -0.4618774 -0.18038532
## Lepilemuridae -0.27611776 0.25220679 -0.7954933 0.20983495
## Lorisidae -0.27017606 0.21479960 -0.7126782 0.13654730
## Pitheciidae -0.31036943 0.07022229 -0.4479488 -0.17578198
## Tarsiidae -0.34849240 0.18409603 -0.7179206 0.01566124
##
## , , scaled_Adult_Length_g
##
## Estimate Est.Error Q2.5 Q97.5
## Aotidae 0.1362992 0.25446620 -0.359433121 0.6720299
## Atelidae 0.4351337 0.07699635 0.269288626 0.5763010
## Cebidae 0.2190944 0.09907545 0.024674468 0.4146046
## Cercopithecidae 0.4316964 0.04052116 0.350499839 0.5105852
## Cheirogaleidae 0.1533433 0.12626433 -0.099440899 0.3971016
## Daubentoniidae 0.2786565 0.22827554 -0.162031344 0.7450809
## Galagidae 0.1389521 0.13449936 -0.125071518 0.4107979
## Hominidae 1.4745848 0.09183447 1.297929357 1.6686765
## Hylobatidae 0.3242058 0.06836725 0.181681223 0.4547250
## Indriidae 0.3060953 0.08397554 0.139098999 0.4755495
## Lemuridae 0.1752057 0.09221433 -0.009475324 0.3546204
## Lepilemuridae 0.2345731 0.25500899 -0.269735495 0.7344529
## Lorisidae 0.2376647 0.20691881 -0.171223875 0.6409735
## Pitheciidae 0.1842543 0.09556931 -0.009304264 0.3681418
## Tarsiidae 0.1347348 0.13859287 -0.134224373 0.4100331
intercepts <- coef(fit.bones)$MSW05_Family[1] # Specifying the first column only
slopes <- coef(fit.bones)$MSW05_Family[2] # Specifying the second column only
summary(fit.bones)$MSW05_Family
## Warning: There were 1 divergent transitions after warmup. Increasing adapt_delta
## above 0.8 may help. See http://mc-stan.org/misc/warnings.html#divergent-
## transitions-after-warmup
## NULL
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Family))+
geom_point(aes(color=MSW05_Family))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Family)) + # slopes for different families
stat_smooth(aes(color=MSW05_Family, fill=MSW05_Family), method="lm", size=0.75,alpha = 0.2)
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
#+ # average slope with SE
coord_cartesian(xlim = c(-2,2), ylim = c(-0.75,1))
## <ggproto object: Class CoordCartesian, Coord, gg>
## aspect: function
## backtransform_range: function
## clip: on
## default: FALSE
## distance: function
## expand: TRUE
## is_free: function
## is_linear: function
## labels: function
## limits: list
## modify_scales: function
## range: function
## render_axis_h: function
## render_axis_v: function
## render_bg: function
## render_fg: function
## setup_data: function
## setup_layout: function
## setup_panel_guides: function
## setup_panel_params: function
## setup_params: function
## train_panel_guides: function
## transform: function
## super: <ggproto object: Class CoordCartesian, Coord, gg>
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Family))+
geom_point(aes(color=MSW05_Family))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Family)) + # slopes for different families
stat_smooth(aes(color=MSW05_Family, fill=MSW05_Family), method="lm", size=0.75,alpha = 0.2) +
coord_cartesian(xlim = c(100,900), ylim = c(-5,20000))
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Family))+
geom_point(aes(color=MSW05_Family))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Family)) + # slopes for different families
stat_smooth(aes(color=MSW05_Family, fill=MSW05_Family), method="lm", size=0.75,alpha = 0.1) +
coord_cartesian(xlim = c(100,500), ylim = c(-5,5000))
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Genus))+
geom_point(aes(color=MSW05_Genus))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Genus)) + # slopes for different families
stat_smooth(aes(color=MSW05_Genus, fill=MSW05_Genus), method="lm", size=0.75,alpha = 0.2) + theme(legend.position="none")
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Genus))+
geom_point(aes(color=MSW05_Genus))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Genus)) + # slopes for different families
stat_smooth(aes(color=MSW05_Genus, fill=MSW05_Genus), method="lm", size=0.75,alpha = 0.2) +
coord_cartesian(xlim = c(100,900), ylim = c(-5,20000)) + theme(legend.position="none")
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Genus))+
geom_point(aes(color=MSW05_Genus))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Genus)) + # slopes for different families
stat_smooth(aes(color=MSW05_Genus, fill=MSW05_Genus), method="lm", size=0.75,alpha = 0.1) +
coord_cartesian(xlim = c(100,500), ylim = c(-5,5000)) + theme(legend.position="none")
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
ggplot(df.primate, aes(x=X13.1_AdultHeadBodyLen_mm, y=X5.1_AdultBodyMass_g, group=MSW05_Genus))+
geom_point(aes(color=MSW05_Genus))+
stat_smooth(method="lm", se=FALSE, size=0.75, aes(color=MSW05_Genus)) + # slopes for different families
stat_smooth(aes(color=MSW05_Genus, fill=MSW05_Genus), method="lm", size=0.75,alpha = 0.1) +
coord_cartesian(xlim = c(100,400), ylim = c(-5,1000)) + theme(legend.position="none")
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
#Conclusion
In this project, we evaluated the the relationship of mass and head-body length through multilevel regression, considering Family and Genus for Primates.
We can see that more often than not there is a positive relationship between mass and head-body lengths within groups such as their Family and Genus. We also see how the lack of data makes the model very uncertain, relative to the scale at which each group exists. For example, small Primates have a much smaller error, but relative to their size, it expresses just as much of an issue.
What may help, is extra information like age, sex, location and captive/non-captive, other body measurments such as forearm-lengths, ontop of more data for each species (not the mean). There can be even more variation within a species, and this makes it difficult to say that this method is reliable. Even though this is adult data, the generalized approach must reduce the complexity of the information.
Nonetheless, this multilevel regression model exposed both the general relationship between the taxonomic groups, which can at least narrow the physical thresholds of body-mass in certain groups.